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Abstract 


Emergent Dirac fermion states underlie many intriguing properties of graphene, and the 
search for them constitute one strong motivation to explore two-dimensional (2D) allotropes 
of other elements. Phosphorene, the ultrathin layers of black phosphorous, has been a subject 
of intense investigations recently, and it was found that other group-Va elements could also 
form 2D layers with similar puckered lattice structure. Here, by a close examination of their 
electronic band structure evolution, we discover two types of Dirac fermion states emerging 
in the low-energy spectrum. One pair of (type-I) Dirac points is sitting on high-symmetry 
lines, while two pairs of (type-II) Dirac points are located at generic fc-points, with different 
anisotropic dispersions determined by the reduced symmetries at their locations. Such fully- 
unpinned (type-II) 2D Dirac points are discovered for the first time. In the absence of spin- 
orbit coupling, we find that each Dirac node is protected by the sublattice symmetry from 
gap opening, which is in turn ensured by any one of three point group symmetries. The spin- 
orbit coupling generally gaps the Dirac nodes, and for the type-I case, this drives the system 
into a quantum spin Hall insulator phase. We suggest possible ways to realize the unpinned 
Dirac points in strained phosphorene. 

Introduction 

Recent years have witnessed a surge of research interest in the study of Dirac fermions in con¬ 
densed matter systems, ranging from graphene and topological insulator surfaces in two-dimensions 
(2D) to Dirac and Weyl semimetals in 3D, which possess many intriguing physical properties 
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owing to their relativistic dispersion and chiral nature. Especially, 2D Dirac fermion states have 
been extensively discussed in honeycomb lattices, commonly shared by group-IVa elements with 
graphene as the most prominent example, for which Dirac points are pinned at the two inequiv¬ 
alent high-symmetry points K and K' of the hexagonal Brillouin zone (BZ), around which the 
dispersion is linear and isotropic. Later on, 2D Dirac points on high-symmetry lines were also 
predicted in some nanostructured materials, including graphynes ** and rectangular carbon and 
boron allotropes. However, the possibility of 2D Dirac points at generic fc-points has not been 
addressed, and such Dirac point has not been found so far. 

Meanwhile, the exploration of 2D materials built of group-Va elements (P, As, Sb, and Bi) 
has just started. Single- and few-layer black phosphorous, known as phosphorene, have been suc¬ 
cessfully fabricated, and was shown to be semiconducting with a thickness-dependent bandgap and 
a good mobility up to ~10^cm^/(V-s), generating intense interest. While 2D allotropes with 
different lattice structures have been predicted and analyzed for the other group-Va elements, 
we note that the puckered lattice structure similar to phosphorene has been demonstrated exper¬ 
imentally for Sb and Bi (down to single-layer) grown on suitable substrates, and been 
predicted for As as well. Motivated by these previous experimental and theoretical works, and 
in view of the ubiquitous presence of the Dirac fermions and the associated interesting physics, 
one may wonder: Is it possible to have Dirac fermion states hosted in such 2D puckered lattices? 
A simple consideration shows that here any possible Dirac point cannot occur at high-symmetry 
points. The reason is that each Dirac point at k must have a time reversal (TR) partner at —k with 
opposite chirality, whereas the BZ of the puckered lattice has a rectangular shape, of which all 
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the high-symmetry points are invariant under TR. Therefore, if Dirac states indeed exist in such 
systems, they must be of a type distinct from those in graphene. 

In this work, we address the above question by investigating the electronic structures of 
group-Va 2D puckered lattices. We find that Dirac fermion states not only exist, but in fact oc¬ 
cur with two different types: one type (referred to as type-I) of (two) Dirac points are located on 
high-symmetry lines; while the other type (referred to as type-II) of (four) Dirac points are lo¬ 
cated at generic fc-points. Depending on their reduced symmetries, dispersions around these points 
exhibit different anisotropic behaviors. Points of each type can generate or annihilate in pairs of 
opposite chiralities, accompanying topological phase transitions from a band insulator to a 2D 
Dirac semimetal phase, and since they are not fixed at high symmetry points, their locations can 
be moved around in the BZ. Particularly, to our best knowledge, the novel fully-unpinned (type-II) 
2D Dirac points are discovered here for the first time. In the absence of spin-orbit coupling (SOC), 
each Dirac node is protected from gap opening by a sublattice (chiral) symmetry, which can in 
turn be ensured by any one of three point group symmetries. The inclusion of SOC could gap the 
Dirac nodes, and in the case of type-I nodes, it transforms the system into a quantum spin Hall 
(QSH) insulator phase. All these properties make the system distinct from graphene and other 2D 
materials. We further suggest that the novel unpinned Dirac points can be experimentally realized 
by the strain engineering of phosphorene. Our discovery therefore greatly advances our funda¬ 
mental understanding of 2D Dirac points, and it also suggests a promising platform for exploring 
interesting effects with novel types of Dirac fermions. 
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Results 


A group-Va pnictogen atom typically forms three covalent bonds with its neighbors. As shown in 
Fig.l for a single-layer phosphorene structure, the P atoms have strong -hybridization character 
hence the three P-P bonds are more close to a tetrahedral configuration. This results in two atomic 
planes (marked with red and blue colors) having a vertical separation of h comparable to the bond 
length. In each atomic plane, the bonding between atoms forms zigzag chains along p-direction. 
The unit cell has a 4-atom basis, whieh we label as Au, Bu, Al, and Bl (see Fig. 1(c)), where U 
and L refer to the upper- and lower-plane respectively. The structure has a nonsymmorphic T>2fe(7) 
spaee group which includes the following elements that will be important in our discussion: an 
inversion center i, a vertical mirror plane perpendicular to ij, and two 2-fold rotational axes C 2 y 
and C 2 z- Note that due to the puckering of the layer, the mirror planes perpendicular to x and 2 
are broken. With the same valenee electron configuration. As, Sb, and Bi possess allotropes with 
similar puckered lattice structures. 

To study the electronic properties, we performed first-principles calculations based on the 
density functional theory (DFT). The details are deseribed in the Methods. The ealeulated geo¬ 
metric parameters of group-Va 2D puckered lattices with D 2 h{ 7 ) symmetry are summarized in the 
Supplementary Information. The obtained struetures agree with the experiments and other theo- 
retieal calculations, 22,29,33 lattice constants a > b, reflecting that the inter-ehain coupling 
is weaker than the coupling along the zigzag chains. The angle 62 increases from ~ 70° for P to 
~ 85° for Bi whereas di remains ~ 95°. The inter-plane separation h as well as the bond lengths Ri 
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and Rs increase by almost 1 A; from P to Bi, while R 2 , the distance between sites of neighboring 
zig-zag chains, increases only slightly, implying that the inter-chain coupling becomes relatively 
more important with increasing atomic number. 

We first examine their corresponding band structures without SOC, whose effect will be 
discussed later. The results are shown in Fig.2. The puckered lattice of P is a semiconductor with 
a bandgap around T-point. From P to Bi, the direct bandgap at F-point keeps decreasing, and a 
drastic change occurs from Sb to Bi where linear band crossings can be clearly spotted along the 
r-X 2 line. Examination of the band dispersion around the two points (labeled as D and D' in 
Fig. 1(d)) shows that they are indeed Dirac points (see Fig.3(a)). Furthermore, along T-Xi line, 
there gradually appear two sharp local band extremum points for both conduction and valence 
bands, where the local gap decreases from P to Sb with the two bands almost touching for Sb, 
yet the trend breaks for Bi. Remarkably, close examination reveals that for Sb and Bi, close to 
each extremum point there are actually two Dirac points on the two sides of the F-Xi line (see 
Figs. 1(d) and 3(b)). The energy dispersions around these Dirac points are shown in Fig. 1(e) and 
1(f), clearly demonstrating the Dirac cone characters. Therefore, two types of Dirac points with 
distinct symmetry characters exist in this system: one pair of type-I Dirac points (D and D') sitting 
on high-symmetry lines and two pairs of type-II Dirac points (near F and F') at generic fc-points. 

The band evolution around F-point from Sb to Bi and the appearance of type-I Dirac points 
in Fig.2 are reminiscent of a band-inversion process. Indeed, by checking the parity eigenvalues at 
F, one confirms that the band order is reversed for Bi around F-point (see Supplementary Informa- 
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tion). For a better understanding, we construct a tight-binding model trying to capture the physics 
around F-point. Since the low-energy bands there are dominated with p^^-orbital character (Fig.2), 
we take one orbital per site, and include couplings along Ri and R 2 in the same atomic plane (with 
amplitudes R and ^2 respectively) as well as nearest-neighbor inter-plane hopping along R^ (with 
amplitude t±) (see Supplementary Information). Written in the basis of {Au, Al, Bu, B^), the 
Hamiltonian takes the form: 


n{k) = 


0 Q{k) 


( 1 ) 


Q\k) 0 

where Q{k) is a 2x2 matrix of the Fourier transformed hopping terms (see Supplementary Infor¬ 
mation). The Hamiltonian (1) can be diagonalized and possible band crossings can be probed by 
searching for the zero-energy modes, which exist when the condition A = t±/[2{ti + ^ 2 )] < 1 is 
satisfied, with two band touching points at (0, ±A:d) where ko = (2/6) arccos(A). The direct gap 
at r can be obtained as A = 2 [tx — 2{ti + 62 )]- Hence this simple model indeed captures the emer¬ 
gence of two Dirac points D and D', along with a transition as parameter A varies: when A > 1, 
the system is a band insulator; when A < 1, it is a 2D Dirac semimetal. The transition occurs at the 
critical value Ac = 1 when the conduction and valence bands touch at F-point and the band order 
starts to be inverted. This corresponds to a quantum (and topological) phase transition, during 
which there is no symmetry change of the system. 


Model (1) captures the trend observed in DFT results. The overlap between-p^ orbitals is 
larger along the R^ bond, hence one expects that t± > ti > t 2 . By fitting the DFT bands around 
F-point, one finds that from P to Bi, t± decreases a lot, while 62 increases and becomes rela- 
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lively more important (Supplementary Information). The result shows that changes 

from (2.50, 0.77,0.33) for Sb to (1.86,0.63, 0.35) for Bi (units in eV). Hence A crosses the critical 
value from Sb to Bi, indicating the band inversion at T and the appearance of two Dirac points. 

The emergence of low-energy relativistic chiral modes is the most remarkable property of 
Dirac points. To explicitly demonstrate this, we expand Hamiltonian (1) around each Dirac 
point, which leads to the low-energy Hamiltonian 

Hr{q) = v^q^ay -\- TVyQya^, ( 2 ) 

where q is the wave-vector measured from each Dirac point, r = ±1 for D and D', cr/s are Pauli 
matrices for the sub-space spanned by the two eigenstates at the Dirac point (apart from the Bloch 
phase factor): |mi) = (0, 0 , 1, —l)/\/2 and \u 2 ) = (1, —1, 0 , 0)/-s/2, and f 2 )/(fi+f 2 ) 

and Vy = by^4:{ti + ^ 2 )^ — are the two Fermi velocities. The form of (2) may also be argued 
solely from symmetry. Compared with graphene, these type-I points are unpinned from the high- 
symmetry points. They can be shifted along r-X 2 (and even pair-annihilated) by varying system 
parameters such as A, although they cannot go off the line as constrained by the symmetries. 
In addition, different from graphene, ^ the dispersion here is anisotropic, characterized by two 
different Fermi velocities. 

Next, we turn to the fully-unpinned type-II Dirac points. The four type-II Dirac points start 
to appear for Sb in our DFT result, located close to the T-Xi line. They can be more clearly seen 
for Bi (see Fig.3(b)). Again the band evolution implies a local band inversion near F and F'. Here 
F and F' (on F-Ai) are the mid-points of the lines connecting each pair of the type-II points. The 
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low-energy bands there are mainly of Pa;-orbital character. To reproduce the fine features using 
a tight-binding model would require more hopping terms. Instead, we construct a low-energy 
effective Hamiltonian around point F (F') based on symmetry analysis. There the Hamiltonian 
is constrained by cr^ which maps inside each pair (labeled by /r = ±1 for F and F'), and by i, 
C2y, C2z, and TR that map between the two pairs. Expansion to leading order in each wave-vector 
component g, gives (see Supplementary Information) 

Hy{q) = wq^ay + (-mo + + miql)a^, (3) 

where q is measured from F (or F'), w, w', mo, and mi are expansion coefficients. Two Dirac 
points appear at (0, ±go) with go = when sgn{mo/mi) = 1, corresponding to a local 

band inversion around q = 0. Further expansion of the Hamiltonian around the Dirac point (0, vqo) 
{p = ±.1) leads to 

Hy = wq^ay + [ 2 pqo{qy - uqo) + tJ.w'q^]a^. (4) 

This demonstrates that the two points at z/ = ±1 are of opposite chirality, as required by cr„. The 
dispersion is highly anisotropic (at leading order, characterized by three parameters: w,qo, and w') 
because the Dirac point is at a generic fc-point with less symmetry constraint, as compared with 
type-I Dirac points. 

Unlike in 3D systems, Dirac nodes in 2D have a codimension of two hence are generally not 
protected from gap-opening. In the absence of SOC, however, the Dirac nodes here are stable 
due to the protection by sublattice (chiral) symmetry between {Ai} and {Bi} (i = U,L) sites, 
which allows the definition of a winding number (i.e., quantized Berry phase in units of tt) 
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along a closed loop £ encircling each Dirac point: = §^Ak- dk/ir = ±.1, where Ak is the Berry 

connection of the occupied valence bands. And for a 2D Dirac point, the sign of (or the ±7r 
Berry phase) is also referred to as the chirality. ® Using DFT results, we numerically calculate the 
Berry phase for each Dirac point and indeed confirm that they are quantized as ±7r. The signs are 
indicated in Fig. 1(d). 

More interestingly, in the puckered lattice with a four-atom basis in a non-coplanar geometry, 
the sublattice symmetry can be ensured by any one of three independent point group symmetries: i, 
C 2 y, and C 2 z- The resulting protection of Dirac nodes can be explicitly demonstrated in low-energy 
models. For example, consider the type-I points described by Eq.(2). There the representations of 
i, C 2 y, and C 2 z (denoted by V, TZy, and Tlz, respectively) are the same, viz., Then the symmetry 
requirement TZyHriqx, = Hri—qx, qy) by C 2 y directly forbids the presence of a mass term 

moz. Meanwhile, since i and C 2 z maps one valley to the other, they protect the Dirac nodes when 
combined with TR {or Uy if it is unbroken), e.g., considering the combined symmetry of C 2 z and 
TR (with representation T = A the complex conjugation operator): {7lzT)Hr{q){7lz'T)~^ — 
Hr{q), which again forbids a mass generation. The underlying reason i, C 2 y, and C 2 z each protects 
the Dirac node is that they each maps between the two sublattices hence ensures the sublattice 
(chiral) symmetry. In comparison, the mi rror plane Oy maps inside each sublattice, hence it alone 
cannot provide such protection. This reasoning is general and applies to type-II points as well. (In 
model (3), i,C 2 y, and C 2 z have representations as by construction, and when combined with TR, 
again each forbids the generation of a mass term ~ See Supplementary Information.) We 

stress that the three symmetries i, C 2 y, and C 2 z each protects the Dirac points independent of the 
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other two. For example, we could disturb the system as in Fig.4 such that only one of the three 
symmetries survives. The corresponding DFT results confirm that the Dirac nodes still exist. Thus 
the crystalline symmetries actually offer multiple protections for the Dirac nodes in the current 
system. 

SOC could break the sublattice symmetry. Hence when SOC is included, the Dirac nodes 
would generally be gapped. For type-I points, treating SOC as a perturbation, its leading-order 
symmetry-allowed form is i/soc = where is Pauli matrix for real spin. This is similar 

to the intrinsic SOC term in graphene, which opens a gap of 2|A| at the Dirac points. For 
the type-II points, we obtain i^soc = VQy'^z^z in model (3) hence a gap of 2go|??| is also opened 
at these Dirac points. Gap opening by SOC is closely related to the QSH insulator phase. 

Here the band topology can be directly deduced from the parity analysis at the four TR invariant 
momenta. This means that only the band inversion at F between the two type-I points contributes 
to a nontrivial Z 2 invariant; whereas that associated with type-II points does not. It follows that 
Sb is topologically trivial since it has only type-II Dirac points, while Bi is nontrivial since it has 
additional type-I points. These results are in agreement with previous studies. 

Breaking all three symmetries i, C 2 y and C 2 z can also generate a trivial gap term at the 
Dirac points, which competes with the SOC gap. For example, this happens when each atomic- 
plane forms additional buckling structure. Nevertheless, as long as the trivial mass term does not 
close the SOC-induced gap, by adiabatic continuation the band topology will not change. 
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Discussion 


Due to their different locations and the associated symmetries, the two types of Dirac points 
here exhibit properties distinct from that of graphene. With preserved sublattice symmetry and 
in the absence of SOC, the Dirac nodes are topologically stable—they can only disappear by pair- 
annihilation between opposite chiralities. This is unlikely for graphene since the Dirac points there 
are pinned at the high-symmetry points. In contrast, the two types of Dirac points here are less 
constrained. Pair-annihilation (pair-generation) indeed occurs during the quantum phase transition 
as observed from the band evolution. 

It is noted that similar type-I points were also predicted in a few nanostructured mate¬ 
rials. Meanwhile the type-II points discovered here are completely new. They are fully- 
unpinned and has highly anisotropic dispersions. With this discovery, now we can have an almost 
complete picture: 2D Dirac points can occur at high-symmetry points, along high-symmetry lines, 
and also at generic fc-points. 

It is possible to have Dirac points, originally sitting at high-symmetry points, to become 
unpinned when crystalline symmetry is reduced due to structural distortions. However, we stress 
that the type-II points here are distinct in that they are realized in a native crystalline structure 
with relatively high symmetry. Only in such a case, we can have a sharp contrast between generic 
fc-points where the group of wave vectors is trivial and the high-symmetry fc-points where the 
group is non-trivial, and accordingly the type-II point can move around (hence fully-unpinned) 
without any symmetry-breaking. More importantly, it is just because that type-II points occur in 
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a state with high symmetry that the Dirac nodes can be protected (in the absence of SOC): as we 
discussed, the various crystalline symmetries ensure the protection of the Dirac nodes from gap 
opening. 

It is remarkable that the two different types of Dirac points can coexist in the same 2D 
material. We emphasize that it is a result of the lattice structure and the valence character of 
the pnictogen elements. Our DFT result indeed shows that even starting from the P lattice, the 
two types of Dirac points can be separately tuned to appear or disappear by lattice deformations. 
For example, we find that the type-II Dirac points can be generated in phosphorene by applying 
uniaxial tensile strains along the ^-direction. The DFT result in Fig.5 and Fig.6 indeed shows the 
band inversion on F-Xi and the formation of four type-II Dirac points around a strain of 16%. 
Since phosphorene has excellent mechanical properties and a critical strains > 25% has been 
predicted, it is promising that the novel strain-induced topological phase transitions and the 
appearance of type-II Dirac points can be directly observed in strained phosphorene. The lattice 
deformation that produces type-I points is discussed in the Supplementary Information. Similar 
scenarios occur for other group-Va elements as well. 

So far, single-layer As, Sb, and Bi in their free-standing form have not been realized yet. Nev¬ 
ertheless, in view of the rapid progress in experimental techniques, we expect that these materials 
could be fabricated in the near future. Especially, for Sb and Bi, the puckered structures have been 
demonstrated by PVD growth down to single-layer thickness on suitable substrates. Besides 
the topological properties, the presence of Dirac states is expected to endow these 2D materials 
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with many intriguing properties for applications, such as the very high mobility, the half-quantized 
quantum Hall effect, the universal optical absorption and etc. Due to the highly anisotropic 
dispersions of these new Dirac points, the electronic transport properties such as the conductivities 
would show strong direction dependence. In addition, since there is no symmetry connection be¬ 
tween the two types of Dirac points, when they are both present, it is possible to independently shift 
each type of points relative to the Fermi level, e.g. by strain engineering, leading to self-doping 
and even the interesting scenario with both electron-like and hole-like Dirac fermions in the same 
system. With the multiple Dirac points with different chiralities, it is possible to further control the 
carriers near different Dirac points for valleytronic applications. 

In conclusion, based on first-principles calculations of 2D allotropes of group-Va elements 
with puckered lattice structure, we predict the coexistence of two different types of Dirac points: 
Dirac points on high symmetry lines and at generic fc-points. In particular, the 2D Dirac points at 
generic fc-points are fully-unpinned, have highly anisotropic dispersions, and are discovered here 
for the first time. Combined with low-energy effective modeling, we unveil the low-energy prop¬ 
erties of these Dirac points. We show that their appearance is associated with the band inversion 
process corresponding to a topological phase transition. The topology/symmetry protection of the 
Dirac nodes is analyzed in detail. Interestingly, because of the unique lattice structure, there is 
a triple-protection of the nodes by three independent point group symmetries. This also implies 
versatile methods to control the locations as well as the dispersions around the Dirac points. When 
SOC is strong, the Dirac nodes are gapped, and in the case of type-I points (such as for Bi), this 
drives the system into a QSH insulator phase. We further show that the topological phase transition 
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and the novel unpinned Dirac points can be realized in strained phosphorene. Our work represents 
a significant conceptual advance in our fundamental understanding of 2D Dirac points. The result 
also suggests a new platform to explore novel types of 2D Dirac fermions both for their fascinating 
fundamental properties and for their promising electronic and valleytronic applications. 

Methods 

First-principles calculations. Our first-principles calculations are based on the density functional 
theory (DFT) implemented in the Vienna ab-initio simulation package. The projector augmented 
wave pseudopotential method is employed to model ionic potentials. Kinetic energy cutoff is set 
to 400 eV and fc-point sampling on the rectangular BZ is with a mesh size 20x20. The minimum 
vacuum layer thickness is greater than 20 A; which is large enough to avoid artificial interactions 
with system images. The structure optimization process is performed including SOC with the 
local density approximation (LDA) for the exchange-correlation energy"^^ and with van der Waals 
corrections in the Grimme implementation. The force convergence criteria is set to be 0.01 eV/A. 
Hybrid functional (HSEOb)"^’ is used for the band structure calculations. 
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Figure 1 (a,b) Top- and side-view of the 2D puckered iattice structure. The green 

shaded region marks the unit cell, (c) Four-atom basis sites in a unit cell, (d) 2D Brillouin 
zone with high-symmetry points. The locations of Dirac points are schematically marked 
by the red dots: two type-l Dirac points at D and D' on r-X 2 ', and four type-ll Dirac points 
around F-Xi forming two mirror image pairs. F {F') on T-Xi iabeis the mid-point of each 
pair. +{-) indicates the chirality of each point, (e) and (f) show the schematic energy 
dispersions around point D and point F respectiveiy, corresponding to the resuit in Fig.3. 

Figure 2 Band structures of group-Va elements with 2D phosphorene iattice structures 
in the absence of SOC. The size of red (green) dots denotes the weight of projection onto 
Pz (Px) atomic orbitais. 

Figure 3 Dispersion around (a) type-l and (b) type-ll Dirac points for 2D puckered Bi. 
Black dots mark the Dirac point locations. Energy is in unit of eV. 

Figure 4 Perturbations which ieave oniy one of the three symmetries preserved: (a) C 22 
survives when appiying E field in ; 2 -direction; (b) C 2 y survives when shifting the two basis 
sites aiong y in the same direction; (c) i survives when shifting the two sites along y but 
in opposite directions. For each case, the Dirac nodes are stiii protected in the absence 
of SOC. The iower panei in each sub-figure schematicaiiy indicate the movement of the 
Dirac points under each perturbation. For the case in (a), the type-ll points do not show 
appreciabie change in their iocations under a weak electric fieid, while the two type-l 
points annihiiate with each other when they meet at F-point. 
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Figure 5 Band structure of phosphorene with uniaxial tensile strain applied along y- 
direction. Band inversion along r-Xi line occurs around a value of 16%. 


Figure 6 Band structure of phosphorene with uniaxial tensile strain applied along y- 
direction. Here the dispersion is plotted along a path crossing point F (as indicated in 
Fig.5) and perpendicular to the T-Xi direction. Two type-ll Dirac points near F can be 
clearly observed at strain value > 16%. 
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